NUMERICAL HOMOTOPIES TO COMPUTE GENERIC POINTS ON 
POSITIVE DIMENSIONAL ALGEBRAIC SETS 
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ANDREW J. SOMMESE AND JAN VERSCHELDE 

Abstract. Many applications modeled by polynomial systems have positive dimensional 
solution components (e.g., the path synthesis problems for four-bar mechanisms) that are 
challenging to compute numerically by homotopy continuation methods. A procedure of A. 
Sommese and C. Wampler consists in slicing the components with linear subspaces in general 
position to obtain generic points of the components as the isolated solutions of an auxiliary 
system. Since this requires the solution of a number of larger overdetermined systems, 
the procedure is computationally expensive and also wasteful because many solution paths 
diverge. In this article an embedding of the original polynomial system is presented, which 
leads to a sequence of homotopies, with solution paths leading to generic points of all 
components as the isolated solutions of an auxiliary system. The new procedure significantly 
reduces the number of paths to solutions that need to be followed. This approach has been 
implemented and applied to various polynomial systems, such as the cyclic n-roots problem. 
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1. Introduction 



Let 



fix) 



fi{xi, 



fN^Xi, 



Xr, 



I Xr. 



denote a system of A^ polynomials on C". Positive dimensional components of the solution 
set of f{x) = are a common occurrence, even when N = n. Sometimes they are an 
unpleasant side show that happens with a system generated using a model, for which only 
the isolated nonsingular solutions are of interest; and sometimes, the positive dimensional 
solution components are of primary interest. In either case, dealing with positive dimensional 
components, is usually computationally difficult. In |^, Sommese and Wampler presented 
a numerical algorithm, which uses auxiliary systems to numerically find sets of solutions of 
the original system. These sets of solutions, which let one numerically decide the dimension 
of the zero set of the original system, include at least the isolated solutions of the original 
system, plus "generic points" of each positive dimensional irreducible component of the 
solutions of the original solutions. Generic points are the basic numerical data which we are 
using to investigate the positive dimensional solution components. 

The algorithm from |3^, which is based on slicing with general linear spaces of different 
dimensions, leads to n auxiliary systems, which must be dealt with. In this paper we present 
an embedding of the system f{x) = into a family of system of polynomials depending on 
2n variables (x, z) G C^", and a large space of parameters. We then single out n + 1 of the 
systems, Si{x, z) for i from n to (here £q{x, z) = is equivalent to the system /(x) = 0) 
obtained by choosing particular values of the parameters, plus a homotopy Hi going from 
£i to £i-i. For simplicity we assume that / is not identically zero. The system £„ = has 
isolated nonsingular solutions. We use polynomial continuation [|I], 0, [^, [^ to implement 
the following algorithm 

Algorithm 1. Cascade of homotopies between embedded systems. 

Input: /, n. system with solutions in C" 

Output: {£i, Xi, Zi)^^Q. embeddings with solutions 



£o := /; 

for i from 1 up to n do 

£i := Emhed{£i-i, Zi); 
end for; 

Zn := Solve (£:„); 
for i from n — 1 down to do 

if,+i := tf,+i + (1 - t) f f' ]; 

Xi := limits of solutions of -ffj+i 

as t — s> with Zi = 0; 
Zi := Hi+i{x,Zi 7^ 0,t = 0); 

end for. 



initialize embedding sequence 

slice and embed 

Zi = new added variable 

homotopy sequence starts 

all roots are isolated, nonsingular, with Zn y^ 

countdown of dimensions 

homotopy continuation 

t : 1 -^ to remove ZiJ^i 

on component 

not on component: these solutions 

are isolated and nonsingular 
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The routine 'Embed' will be defined in the next section. In section three we present a 
worked out example of the algorithm. Section four contains the mathematical background 
needed to prove our main results: 

1. if i is the largest integer with A", nonempty, then the dimension of /^^(O) is i; and 

2. given any irreducible component W of /~^(0) of dimension i, then, counting multiplic- 
ities, Xi contains deg(W^) generic points of W. 

The applications described in section five illustrate the performance of the new procedure. 
We end this paper with some directions for future research. 



Acknowledgements. Both authors would like to thank MSRI for providing the stimulating 
environment where our collaboration began. We would also like to thank the Center for 
Applied Mathematics, and the Duncan Chair of the University of Notre Dame for their 
support of our work. The authors are grateful to Goran Bjorck for pointing at the results 
of Uffe Haagerup. The second author is supported by a post-doctoral fellowship at MSRI, 
and supported in part by the NSF under Grant DMS - 9804846 at the Department of 
Mathematics, Michigan State University. 



2. An Embedding of a Polynomial System 

Throughout this paper we work with algebraic functions. This allows the possibility of 
using rational functions, and not just polynomials. This extra flexibility will be needed 
in a sequel, where, even though we start with a system of polynomials on C", it becomes 
necessary to work with rational functions on a Zariski open set of an associated Euclidean 
space. We assume in what follows that locally we have the same number of equations as 
unknowns. A procedure for reducing to this case is presented in [^. 

We refer to [Q for a discussion of generic points. 

Given a system of algebraic functions 



fix) 



fn{x) 



(2) 



on a connected algebraic manifold X of dimension n embedded into a Euclidean space C"^ 
we have the following basic embedding into a family of systems. First we restrict n linear 
functions on C"^ to X. Thus for j from 1 to n we have 

Lj{x) := aj + ttj^iXi -\ h Uj^aXa, (3) 

where Xi is the restriction of the i-th coordinate function of C"^ to X. By abuse of notation 
we let Lj e C^+^ denote 

( «j «i,i • • • HA ) • (4) 
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We fix linear coordinates zi, . . . ,Zn on a. complex Euclidean space C". Then we have the 
system of equations 



£^{f){x,z,Al,... ,Ai,Li, 



■.U) 



fli + ai^ia;i H h ^i.^x^ + Zi 

tti + ttiixi H h aiAXA + Zi 



(5) 



where we let 



A,: = 



Ai,i 



An,; 



(6) 



We often refer to £{f){x, zi, . . . , Zn, Ai, . . . , A„, Li, . . . , L„) by Sn or £„.{,!)■ Further we let 
£i or £i{f) denote £{f){x, zi, . . . , Zi, Ai, ... , Aj, Li, . . . , Li) on €"+*. Note that £o is just /; 
and that the solutions {x,zi, . . . ,Zi) E X x C^ of 

£i{x, zi,... ,Zi,Ai,... , Ai, Li, . . . ,Li) =0 (7) 

are naturally identified with the solutions {x,Zi, . . . ,Zi,0... , 0) G X x C" of the system 



£i{x, Zi, . . . , Zn, Ai, . . . , Aj, 0, . . . ,0, Li, ... , Lj, 0, . . . ,0) 

Zi+l 



We let Y denote the space C"^'^^+-'^^ x C"^" of parameters 
ai ail ■ ■ ■ a-i,A Ai i ■ ■ ■ A^.i 






0. 



O'Ti CLnA ■ ■ ■ O-n 4 Ain " " " ■^n-j 



^ (pnx(A+l) ^ (pnxn 



(9) 



for these systems. We have used the transpose of the Xij for convenience in describing the 
stratification Yq C Yi C ■ ■ ■ G Yn oi the space Y given by defining Y^ := Y; and Yi, for 
z = 0,...,ra — 1, as the subset of Y with the coordinates 



O'i+l Cti+l,! 



^n ^riA 



0'i+l,A ^1,1+1 ■ ■ ■ Xn,i+l 



0'n,A Ai„ 



'^n,n 



(n— «)x(A+l) y, (p(n— j)xn 



(10) 



set equal to 0. Thus using the identification @ above, we can regard Yi as the parameter 
space of the systems of equations 



£i\x, Zi, . . . ,Zi, Ai, . . . , Aj, Li, . . . , Li) — 0. 



:iii 



Hi{x,zi,... ,Zi,t) :-- 
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We consider homotopies Hi, for i from n to 1 and t from 1 to 0, defined by: 

fiix) + E}=i >^i,jZj + tXi^iZi 

JnyXj + Z^j=l ■^n,jZj + tAn,iZi 

Li{x) + zi 

Li-i{x) + Zi_i 
tLi{x) + Zi 

with the convention that if i = 1, we mean 

fi{x) +tAi,iZi 

tLi{x) + Zi 

Thus Hi{x, zi, . . . ,Zi,l) = Si and Hi{x,zi, . . . , Zi, 0) = is 

f Si-i{x,zi,... ,Zi-i) = 
I Zi =0 



Hi{x,zi,t) :-- 



Note that using this convention, Hi can be rewritten as tSi + {1 — t) 

Lemma 2. There is a nonempty Zariski open set U of points 
tti ail ■ ■ ■ o,i,A -^1,1 ■ ■ ■ ^n,l 



Zi 



(12) 



(13) 



(14) 



(15) 



such that for each i = 1, . . . ,n, 

1. the solutions of Si{x, Zi, . . . ,Zi) = with {zi, . . . , Zi) y^ are isolated and nonsingular; 

2. given any irreducible component W of f~^{0) of dimension i, and counting multiplici- 
ties, the isolated solutions of Si{x, Zi, . . . , Zi) = with {zi, . . . ,Zi) = 0, contain deg(W^) 
generic points of W ; and 

3. the solutions of Si{x, zi, . . . , Zi) = with {zi, ... ,Zi)^0 are the same as the solutions 
of £i{x, zi, . . . ,Zi) = with Zi ^ 0. 

Proof. For a fixed i E {0, . . . ,n} consider the system of equations 



f EUaiJ, + r,=i f3i,z, = 

l^j=l'^n,jjj + l^j=lPn,jZj = U 

ii + Y.f=iii,jXj + Y.'j=iSi,jZj = 



(16) 
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By Bertini's theorem (see section four for a convenient form of this result), there is a 
nonempty Zariski open set U of 

(a,/3,7,(5) e C"^" X C"^^ X C^(^+^) x C^^ 



(17) 



where 



a 



«i,i 



"n,! 



ttl,n 



Ctr. 



(3:- 



and; 



7 



7i 7i,i 



7j 7i,i 



7l,A 



7i,A 



/3i,i 



Pn,i 



(19) 



such that, if not empty, the zero set Z of the above system of equations on X x C* minus the 
set of common zeroes of /i, . . . , /„, zi, . . . , Zj is smooth and of dimension 0. Left muhiplying 
this n + i vector of equations with an invertible (n + i) x (n + i) matrix G of the form 



A 
B 



(20) 



where ^ is an n x n matrix and i3 is a i x i matrix, resuhs in the equivalent system of 
equations 



A-a- f + A- p 



i3-7 



1 

Xi 
Xa 



+ B-6 



Z\ 

Zi 
Z\ 



(21) 



0. 



Thus we can assume that the set V is invariant under this action by the matrices G. Thus 
we have a nonempty Zariski open set \J of 



fll ^1,1 ■ ■ ■ 0'\,A -^1,1 ■ ■ ■ -^n,! 



Q-n On,l ■ ■ ■ On, A ^\,r, 



A 



n,n 



p (pnx(A+l) ^ (pnxn 



(22) 



such that for each i the equivalent system is of the form 

f fl{x)+El=lK^Z^ = 



/n(a;) + Ei=l -^n,i^i = 

«i + E/=i oi jXj + 2;i = 



(23) 



(li ~r 2-^j^i ^i,j'^j "T -^^ 
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Zi) 7^ 0. Thus we have the first assertion 

with 



and has smooth nonsingular zeroes when [zi, . 
of the Lemma. 

To see the second assertion of the Lemma, note that the solutions of Sj 



[zi, 



are naturally identified with the solutions of the system 

fiix) = 



fn{x) 



«i + Ej=i a 



«j-^j 










(24) 



The second assertion follows now from the Algorithm in [^, §3.1]. 

We prove the third assertion of the lemma by induction on z. If z = 1, then it is a tautology 
that the solutions of Si{x, Zi, . . . , Zi) = with {zi, . . . ,Zi) 7^ are the same as the solutions 
of Si{x, zi, . . . , Zi) = with Zi 7^ 0. 

So we can assume that the result is true for k < i where i > 1. Note that a solution of 



Si = with Zi = but (zi, 



7^ is a solution of the system 



Si-i{x, zi, 



! Zi-l) 

Liix) 







0. 



(25) 



Since the solutions of Si^i = with Zi_i 7^ are isolated and nonsingular, a generic choice 
of Li{x) will not be zero on any of the solutions. But this means that the solutions of £^j = 
for generic Li{x) will have no solutions with Zi = 0. □ 

Note that if we choose y generically in Y, then we have chosen the associated i/i generically 
in Yi for each i from n to L Thus we can assume that with an initial generic choice of 
parameters y, the behavior for each of the systems £^i = is the behavior we expect with 
a generic choice of parameters on Yi. One minor point remains. Given a generic choice of 
parameters y & Y, generic behavior might not occur for the homotopy Hi with t G (0, 1]. 
This is easily dealt with by a trick of Morgan and Sommese |2^, ^, 3D . 

Assume that we have chosen the parameters y for the systems in the nonempty Zariski 
open set U of 



di (i\ 1 



0>nA 



«1,A Ai 



0'n,A Ai n 



A 



ra,l 



^ra,n 



^ (j-nx(A+l) 



X 



(26) 



We want our homotopies Hi{x,zi, . . . , Zi,t) = to define algebraic sets that are fiat over 
a Zariski open set containing (0, 1]. Numerically this means that we want to have generic 
behavior for each t G (0, 1], i.e., the same number of isolated points, no components that do 
not correspond to components in any other fiber. Exploiting generic flatness in this way is 
the underlying approach of the work of Morgan and Sommese [E3, 123, RD| , to which we refer 
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for more details. Since for all but a finite number of values of t G C generic behavior occurs, 
we can conclude that given a set of parameters y & U as above, the homotopy 

f I ^\ I v^i— 1 

h 



^1^,1X^1 ^1) 



1 Zi,t) .- 



{x,zi, 



vt) 



fiix) + 



E}=i >^i,jZj + trjXi^iZi 



Li{x) + zi 



(27) 



Li_i{x) + Zi-i 
triLi{x) + Zi 

defines an algebraic set flat over a Zariski open set of C containing (0, 1], for all but a finite 
set of r] G C with \t]\ = 1. Using openness of flatness, we can absorb the t] into the parameters 
we use. We still have a dense open set of general parameters, but the set is only Zariski open 
in the underlying real algebraic structure. 

By Lemma ^ if the Li, Aj are chosen randomly, then with probability one, Sn{x, z) = has 
a solution with z = only if f{x) is identically zero on C". 

Recall from the introduction that for i from n to 1, 

1. Zi denotes the solutions to Si = with Zi ^ 0; and 

2. Xi_i denote the limits with Zi^i = of the paths of the homotopy Hi{t), from t = 1 to 
t = 0, starting at points of Zi. By convention, the condition Zi^i = is empty when 
i = 1. 

Theorem 3. Let f be as above. Assume that f is not identically zero and that the Aj are 
chosen generically. Ifi is the largest integer with Xi nonempty, then the dimension o//~^(0) 
is i. Moreover given any irreducible component W of f^^{0) of dimension i, then, counting 
multiplicities, Xi contains deg(W) generic points ofW. 

Thus our algorithm achieves the same numerical goal of the algorithm of [^ , but much 
more efficiently. 

Proof of Theorem |^. We use the notation of the proof of Lemma |[ By Lemma H, it follows 
that there is a nonempty Zariski open set of points y & Y, such that for each i from n to 
1, all elements of the set Zi of solutions of Si{x,zi, . . . , Zi) = with Zi ^ are nonsingular 
and isolated. Thus for a dense set (Zariski open in the underlying real algebraic structure of 
7/ G y) the paths over t G (0, 1] are smooth if they start at nonsingular isolated solutions of 
Si = 0. By Lemma |^ and our random choice oi y G Y, each nonsingular isolated solutions 
(x*, zl, . . . , z*_-^^, 0) oi Si^i = must start from a nonsingular isolated solution (x', z[, . . . , z'^) 
of S, = 0. 

If z'^ = for a solution, then {z[, . . . ,z[) = by Lemma |[ Thus f{x') = and hence 
Hi{x', 0, . . . , 0, t) = identically in t, and hence {x*, z^,. . . , 2;*_^, 0) = {x', 0, . . . ,0). Thus 
we know from Lemma ^ that for a Zariski open set of Yi the solutions of the system Si = 
with Zi = (and hence [zi, . . . , Zi) = 0) are on irreducible components of the algebraic 
set f{x) = of dimension i. Thus (x*, zl, . . . , 2;*_^, 0) = (x', 0, . . . ,0) is on an irreducible 
component of the algebraic set of points with Si-i = 0, which has dimension > 1. This 
contradicts {x*,zl, . . . ,2*_i) being an isolated nonsingular solution of Si-i = 0. Thus the 
isolated nonsingular solutions ofSi-i = must be endpoints of paths of homotopy Hi starting 
at points of Zi. 
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Using Lemma ^ the above argument can be modified to show that isolated, but possibly 
singular, solutions oi {x*,zl, . . . , z*_-^^, 0) of £i_i = must start from a nonsingular isolated 
solution {x',z[, . . . , z'^) oi £i = with {z[, . . . , z'^) ^ (and hence by Lemma with z\ ^ 
0). D 

The property of a generic system of the form £n{f) that it has isolated nonsingular solutions 
is useful. In this direction, see fM . 



3. A Worked Out Example 
Consider the polynomial system f{x) = and start system g{x) = 0. 

fi^) = I 2/ 2 f'"'', ^ n 9ix) = I i " ^1 = j; (28) 

^ ^ [ 3^1(3^2 + 2:1) = ^ ^ [ X2 - C2 = ^ ^ 

To solve f{x) = we trace D = 3 x A = 12 solution paths starting at the solutions of 
g{x) = 0. One path diverges to infinity, three paths converge to (0,0), and the remaining 
eight paths end at the solution component xi = 0, for some X2 7^ 0. The condition numbers 
at the end of the paths do not allow us to decide which solution is isolated. 

To embed f{x) =0 we take a random hyperplane L{x) = ao + ^n^i + 0'2a;2 = and choose 
two random complex constants, Ai and A2: 

{xlx2 + XiZ = 
xjixl + xi) + X2Z = (29) 

ao + aiXi + 02X2 + -2 = 

We can solve this system by tracing D = 12 solution paths, using a standard linear homotopy 
with the equations g{x) = 0, z — 1 = as start system, with g{x) = as in (pSf). Five paths 



diverge to infinity. Two paths go to the same solution with 2; = 0, which reveals the degree 
of the solution component xf = 0. Note that geometrically this component corresponds to a 
pair of lines. The five remaining paths go to regular solutions with 2; 7^ 0. 

To compute the possible remaining isolated solutions, we trace five solution paths starting 




at the five regular solutions of the system (|29D . In going with t from 1 to 0, we use the 
homotopy 

/ f xlx2 + XlZ = 

Hi{x,z,t) =t\ < xl{xl + Xi) + X2Z = 

\[ ao + aiXi + 02X2 + z = J \[ z = J (30) 

Three paths converge to (0, 0), two paths go to solutions on the component with xi = and 
X2 7^ 0. End games are still needed to decide whether the solutions are isolated. 

With our new method, 17 solution paths instead of 24 were traced, as 24 would have been 
the number of paths with an iteration of the procedure in ||35|| . 

Polyhedral root counting methods provide a generically sharp root count for polynomial 
systems. In particular, the mixed volume of the Newton polytopes of the system equals the 
number of roots in (C*)", C* := C \ {0}, for a system with generic coefficients. When the 
system has only few monomials with nonzero coefficients, then the mixed volume provides a 
much lower root count than the Bezout bounds based on the degrees of the polynomials. 
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For our type of applications, the distinction between solutions with z = and z y^ 
is instrumental in identifying components of solutions. This difference does depend on the 
values of the coefficients of the original system and is not neglected by the ordinary mixed 
volume. So we will not miss any solutions with z = 0, but we may miss solution components 
for which some Xi = 0. Fortunately, extensions of the polyhedral methods that allow to 
count and compute all affine roots (that is in C" instead of (C*)") are covered amply in the 
literature (see ^, 0, ^, |25|, |2|, |3[ and g). The key idea [|5[ is to add a random 



constant to every equation to shift the roots with zero components away from the coordinate 
axes. In removing these constants by continuation, all affine roots lie at the end of some 
path that starts at a root in (C*)"". 

Consider again the polynomial system f{x) = in (^). Because of the first equation 
we immediately see that there cannot be any solution with all components different from 
zero. The direct application of polyhedral methods does not yield anything, since the mixed 
volume for f{x) = equals zero. With affine polyhedral methods, we consider the system 

where 71 and 72 are randomly chosen complex numbers. Here the mixed volume equals 
five. Note the difference with the total degree D = 12. Letting the 7's go to zero, three of 
the five paths converge to the origin, and the other two paths go to other solutions on the 
component. 

To obtain information about the components, with polyhedral methods we consider the 
embedding 

{xlx2 + 7i + Xiz = 
xl{xl + xi) + -f2 + X2Z = (32) 

ao + aiXi + a2X2 + z = 

This system has mixed volume equal to seven. We use this as start system to solve it with 
7i = and 72 = 0, following the paths that start at the seven solution paths of Si{x, z) = 0. 
From the seven paths, five paths go to solutions with z y^ 0, and the other two paths go to 
the component ending with z = 0. 



4. Bertini's Theorem and a Local Extension Theorem 

Here is a weak, but convenient form of Bertini's theorem, e.g., Fulton Jl^ Example 12.1.11]. 
For a further discussion of Bertini theorems, see also |^, §1.7]. 

Theorem 4 (Bertini). Let X be an algebraic manifold of dimension n, e.g., complex Eu- 
clidean space or complex projective space. Let £1, . . . , £„ be line bundles on X . For i from 1 
to n, let {sjjij = 1, . . . , Tj} be a set of sections of Ci. Let Bi denote the set of common zeroes 
of the sections {sij\l < j < Tj}; and let B := Uji?j. Then given general real or complex 
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numbers {Xi,j\j = 1 



, . . . , Tj, ? - 


=1,... 


, n}, all solutions on X - 

r -ri 

}_^Xi,jSij = 
1=1 


- B, of the system of equations 




' 


) ] ^n,jSn,j = 

I i=i 


(33) 



are isolated and nonsingular. 

In the proof of Theorem ^, we need an extension theorem, giving conditions ensuring that 
if we have a system of equations f{x, y) = depending on parameters y with a multiphcity 
fi isolated solution x* of the system for some point y* in the parameter space, then there 
is a neighborhood U of x* such that for points y near y* in the parameter space, there are 
/i solutions counting multiplicities of the system f{x,y) = 0. Special cases are well known 
in the context of all polynomial systems, but we do not know a general reference in the 
numerical analysis literature. This sort of result is standard for algebraic geometers in the 
algebraic context, or within the German school of several complex variables in the complex 
analytic context. So we are simply explaining why this sort of result follows immediately from 
standard results in these fields. We work locally. All open sets are in the usual Euclidean 
topology, i.e., not in the Zariski topology. 

A possibly nonreduced complex analytic space, X, is said to be a local complete intersec- 
tion if given any point x G A", there is a set W C A:", open in the usual complex topology, 
that contains x, and such that 

1. there is an embedding cf) : U ^ B oiU into an open ball in 5 C C^ for some N : = 
n + m > 0; 

2. the ideal of the complex space 0(W) is defined by holomorphic functions gi, ■ ■ ■ ,gn', and 

3. the dimension of the maximal dimensional irreducible component of X through x is m. 

Local complete intersections are very well behaved. One elementary, but important fact 
about them is that all the irreducible component of X through a given x E X are equal 
dimensional. A less elementary, but important result is that for surjective holomorphic 
maps (respectively algebraic morphisms) from possibly nonreduced complex analytic spaces 
(respectively possibly nonreduced algebraic varieties) to complex manifolds (respectively 
algebraic manifolds), flatness of the morphism is equivalent to the openness of the morphism 
(which is equivalent to all the irreducible components of all the fibers of the morphism 
restricted to a given connected component of X being equal dimensional). This last result 
is a standard result for complex analytic spaces X, whose structure sheaf consists of local 
rings which are Cohen-Macaulay, e.g., Fischer []T6|, page 158]. Local complete intersections 
are among the simplest examples after manifolds of complex analytic spaces with Cohen- 
Macaulay local rings. 

Let X be a connected n-dimensional complex manifold. Let F be a connected r/i-dimensional 
complex manifold. Let 

fi{x,y) 



f{x,y) 



fn{x,y) 







(34) 
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be a system of n holomorphic functions. Let x* be an isolated solution of the system 
f{x,y*) = for a fixed value y* G Y, i.e., assume that there is an open set O C X 
containing x* with x* the only solution of f{x,y*) = on O. Assume that the multiplicity 
of X* is /i. This number, which is 1 exactly when the Jacobian of f{x, y*) is invertible at x*, 
is equal to 

duncOx\.*/J{fi{x,y*), . . . Jn{x,y*)) (35) 

where Ox\x* is the local ring of convergent power series on X centered at the point x*, and 
J{fi{x, y*),. . . , fn{x, y*)) is the ideal in Ox\x* generated by the functions fi{x, y*),. . . , fn{x, y* 

Lemma 5 (Local Extension Lemma). Let X , Y , f, x* , and y* be as above. There are open 
neighborhoods U of x* & X and V of y* E Y such that for any y E V there exist /i isolated 
solutions (counting multiplicities) of f{x,y) = on U. 

Proof. Let diin[(^x* ,y*)} Z denote the dimension of an analytic set Z C X x F at the point 
(x*, y*). Let X' denote the zero set of /i, . . . , /„ on X x Y. Since there are n-functions, the 
dimension of each irreducible component of X' is at least m, and in particular 

dim.(^x*,y*) X' > m. (36) 

Since x* is isolated we have that divn^x* ,y''){,X x {y*}) fl X' = 0. Since X xY is smooth, we 
have 

dim(^*„y.)(X X {y*}) r\X' > dim(^..,,y*) X x {y*}) + dim(^.„y.) X' - dim(^*,,j^*) X xY 

= n + dim(a;.^^j^.) X' — n — m. 

Thus we conclude that 

dim(j.. y.) X' = m. (37) 

Thus the union X" of the components of X' passing through (x*, y*) has pure dimension m. 



It follows, e.g., use Gunning I^O], Theorem 16], that there are neighborhoods U oix* E X and 
V oiy* eY such that the projection vr of X"r\ {U x V") to \^ is proper and finite. By shrinking 
[/ and l^ we can assume that X"n([/x\/) =X'n{UxV) and that X'n([/x {y*}) = {x*,y*). 
Now X' n {U X V) is a local complete intersection, and thus since the map iix' : X' ^ V 
is proper with finite fibers, we conclude by the discussion before the discussion of flatness, 
that this map is flat. Thus the direct image S := 7r^,((9x') of the structure sheaf of X' is 
locally free, e.g., Fischer |jl6[ Proposition 3.13]. On U x V, X' is defined by the functions 
fi{x,y) = 0, for z = 1, . . . , n. By definition, at a point {x',y') G {U x V) H X' we have 

Ox'\{{x',y')} is 

0{uxv)\{{x',y')}/J{fi{x,y),... Jn{x,y)) (38) 

where 0(^uxv)\{ix',y')} is the local ring of convergent power series onUxV centered at {{x',y')} 
and J{fi{x, y),... , fn{x, y)) denotes the ideal in 0{jj^v)\{{x',y')} generated by /i, . . . , /„. 

The statement that £ is locally free, is equivalent to ranks of £ at different points of V 
being equal. The rank of £^ at a point y' eV is by definition 

dimc^y := £/ (my ®Oy|^, £) (39) 

where iriy' is the maximal ideal generated Oyy consisting of convergent power series vanishing 
at y' . Thus comparing to the definition of the multiplicity of the solution x* of /(x, y*) being 
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fi, we see that the rank of S is fi. Unwinding the definition of vr*, Sy> for any fixed y' ^ V is 
the direct sum of the modules 

Oxi.'/Jifiix,y'),...,Ux,y')), (40) 

with index set the set of distinct points x' E U with f{x', y') = 0. This proves the assertion 
of the Lemma. □ 

5. Applications and Computational Experiences 

We have conducted a systematic set of experiments with PHCpack [^ on four case studies 
of famihar polynomial systems. Little symbolic manipulation of polynomials is needed to 
set up homotopies derived from the embedding presented in this paper. 

To avoid wasting paper, we have omitted the algebraic formulations of the systems, which 
can either be found electronically in the database maintained on the web sites of the second 
author, or can be consulted in the cited literature. Although the systems are academic 
examples and only interesting for benchmarking purposes, we try to indicate the relevance 
to their apphcation fields. 

Unless stated otherwise, all reported timings concern a 450MHz Intel Pentium machine 
with 1Gb of main memory and 1Gb swap space, running Debian GNU/Linux. We observe 
that the initial stage is always the most expensive one. Dealing with components of solutions 
is a much harder problem than just approximating the isolated solutions. 

With a faster computer one can attack larger problems, but also the quality of the soft- 
ware matters. Concerning this latter aspect we want to point out that no special-purpose 
software has been developed that exploits the structure of the embedded systems. The con- 
tinuation could for instance go faster if one eliminated explicitly some variables using the 
linear hyperplanes that were introduced in the slicing. 

5.1. A planar four-bar mechanism. Four-bar mechanisms are ubiquitous in mechanical 
design. The 4- variable polynomial system that was derived and solved in |^ has total degree 



256 and lowest multi-homogeneous Bezout bound 96. The mixed volume equals 80, avoiding 
the calculation of zero component solutions. There is a solution component of dimension 
two, with sum of degrees equal to two. The lowest multi-homogeneous Bezout bound of £2 
equals 240 whereas the mixed volume equals 96. 

We summarize the computational experiments in Table |l|, whose format goes as follows. 
The cascade of homotopies starts with solving the start system g by polyhedral methods since 
mixed volumes are sharper than the bounds based on the degrees. To solve the system Sk we 
must trace ^paths solution curves. This number is partitioned into ones on the component 
{z = 0), regular finite solutions (z 7^ 0), and diverging ones (— *> 00). Observe the take over 
of the {z 7^ 0)-solutions to the next row. The number in the column with heading z ^ on 
the Sq-tow is the number of isolated solutions of the original system. In the last column we 
list cpu times. 

The information in Table |1| is useful to make some rough comparisons with other homotopy 
methods. For instance, if we are only interested in the isolated roots, the cost would be of 
the same magnitude as on the Sq-tow. To compare with the procedure of Sommese and 
Wampler |^ , we could add up the entries in the columns with headings z = and ^ 00 to 
the ^paths in the next row, simulating the process of no solution recycling. 
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system 



(7, start 
£1 

£, 



total 



^paths 



96 
96 
68 
56 



316 








2 

12 
20 



34 



z^Q 



96 
68 
56 
36 



256 



00 





26 







26 



cpu time 



10s 880ms 

4m 47s 830ms 

18s 370ms 

lis 10ms 



5m 28s 90ms 



Table 1 . The number of paths and timings for a planar four-bar mechanism. 



We must note that there exists an isotropic formulation of this problem [37] for which the 
mixed volume is an exact root count for the 36 isolated solutions. The black-box solver of 
PHC needs only lis 710ms to solve that system. 



5.2. Constructing Runge-Kutta formulas. The following system is due to Butcher |TD 
and was used as test problem in [^ . The system has seven variables and a three-dimensional 
component. Its total degree is 4608 and mixed volume equals 24. After slicing and embedding 
thrice, the resulting 10-variable system has mixed volume 247. We summarize the running 
statistics in Table ^ 



system 


#paths 


z = 


Zy^O 


-^ 00 


cpu time 


g^ start 


247 





247 





17m 13s 770ms 


^3 


247 


3 


193 


51 


9m 57s 60ms 


S2 


193 


15 


161 


17 


Im 59s 470ms 


Si 


161 


11 


72 


78 


4m 57s 600ms 


£0 


72 





4 


68 


2m 24s 650ms 


total 


920 


29 


677 


214 


34m 32s 550ms 



Table 2. The number of paths and timings for Butcher's system. 



5.3. Modular equations for special algebraic number fields. The authors of |T^ dis- 
cuss the finding of the solution components of the system that was derived in [|r^. That 
system has a component of dimension two, and in |Tl| the authors reasoned that any three 
of the four polynomials generate the same ideal. The challenge posed in [^ to computer 
algebra was to verify this claim. 



Although we believe that nowadays, more than 10 years after ||TT|, some computer algebra 
packages might be strong enough to compute a Grobner basis for this problem, we would 
like to indicate how one could answer the challenge numerically. The approach we propose 
is to replace any of the four equations by the same random hyperplane. This operation cuts 
the dimension by one. So we have to embed each of the four systems with the same random 
constants. We solve one of them, this takes about 2 minutes, and check whether the generic 
points we find satisfy the other equations. 

5.4. On Fourier transforms: the cyclic n-roots problem. In [^ an example was given 
which stems from the problem of finding all "bi-equimodular" vectors x G C", i.e.: all x with 
coordinates of constant absolute value such that the Fourier transform of x is a vector with 
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coordinates of constant absolute value, see also and 0. This system was popularized 
in |jl2[, and is by far the most notorious benchmark problem in polynomial system solving. 

m and 



Other references are 
In 



14 



26| , the following conjectures of Ralf Froberg are mentioned. If n has a quadratic 



divisor, then there are infinitely many solutions. If the number of solutions is finite, then 
this number equals ( "■_ , ) . Uffe Haagerup |^ has proven that for n prime, the number 



of roots is always finite, and equals indeed L" ^-.,2 . 

In Table ^ we summarize the results of the embedding for the cyclic 8-roots problem, 
which has a solution component of dimension one. The original system has mixed volume 
equal to 2,560. The embedded system has mixed volume 4,176. See Table ^ for the summary 
of the runs. We remark that if one is only interested in the isolated roots, one only has to 
trace 2,560 solution paths. 



system 


T^paths 


z = 


Zy^O 


-^ 00 


cpu time 


g, start 
Si 
So 


4,176 
4,176 
3,975 




144 
495 


4,176 
3,975 

1,152 




57 

2,328 


Ih 30m 26s 800ms 
Ih 10m 3s 930ms 
7h 7m 55s 580ms 


total 


12,327 


639 


9,303 


2,385 


9h 48m 26s 310ms 



Table 3. The number of paths and timings for cyclic 8-roots. 



Before developing the embedding method, we attacked the cyclic 8-roots problem first 
with the slicing method of Sommese and Wampler. The computation of generic points on 
the solution component required the tracing of 10,940 solution paths (which takes now 4,176 
paths). For this task, the black-box solver of PHC needed 4 days and nights 15h 6m 49s 392ms 
cpu time on a SUN workstation. Although the Linux machine we used is somewhat faster 
than that SUN workstation, it certainly does not speed up things 50 times! 

We end this section with another illustration on the importance of having good polynomial 
equations for the same problem. Thanks to a trick of John Canny, described in [|I^, there 
exists a reduced version of the cyclic n-roots problem. This version reduces the number of 
roots by eight and yields a significant saving in the path tracking. For instance, the mixed 
volume of the reduced cyclic 8-roots problem equals 320, instead of 2,560 for the original 
problem. Table H summarizes the computational experiments. 



system 


#paths 


z = 


zj^O 


-^ 00 


cpu time 


g, start 
Si 

So 


775 
775 
744 




18 
445 


775 
744 
144 




13 

155 


10m Is 860ms 
11m 23s 680ms 
14m 52s 340ms 


total 


2,294 


463 


1,663 


168 


36m 17s 860ms 



Table 4. The number o 



paths and timings for reduced cyclic 8-roots. 



Finally, this reduction also allows us to treat the cyclic 9-roots system that has a two- 
dimensional component of solutions. The sum of the degrees of the reduced cyclic 9-roots 
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problem equals two. The mixed volume of embedded system E2 equals 4,044, computed 
in cpu time 2h 48m 4s 320ms. The polyhedral continuation for the start system requires 
Ih 4m 5s 130ms and the path following to solve E2 takes 2h 15m 4s 540ms. In total, the 
black-box solver of PHCpack needs 6h 7m 45s 270ms to solve 82- 

6. Conclusions and Future Directions 

In the paper we have presented an embedding of polynomial systems to compute generic 
points on components of solutions. The homotopies recycle solutions while moving down to 



the isolated solutions. The major advantage compared to the method in |35| is that fewer 
solution paths diverge. We provide practical evidence for the efficiency of the embedding 
technique. 

Besides the reported progress on homotopy methods, we want to point out that our em- 
bedding technique might be of independent interest, in the sense that it allows us to focus 
on one component of a particular dimension. The main novelty of the embedding consists 
in the fact that the solutions of the embedded system that do not lie on the focused com- 
ponent provide meaningful information about the components of lower dimension. We hope 
that besides homotopy methods, other solvers for polynomial systems could benefit from our 
embedding technique. 

It is sometimes said that every scientific problem solved generates three other questions. 
Here we briefly indicate three different directions for future research. First of all, the em- 
bedding leads to higher-dimensional systems with special structure. It would be worthwhile 
to exploit this structure to obtain a more efficient mixed-volume computation to alleviate 
the complexity of the flrst stage. Secondly, the computed generic points on the components 
are an interesting starting point for future computational explorations of the solution com- 
ponents. Note that continuation is a very natural tool to scan the component while wiggling 
the added hyperplanes. Lastly, in the flnal stage of the cascade of embedded homotopies, 
it remains difficult to separate solution paths converging to components from other possi- 
bly ill-conditioned isolated solutions. End games that produce certiflcates in the form of a 
random hyperplane through a nonisolated solution will have to be developed. 
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